%% 参数设置
clear; clc;

a = 0.85;          % QO-IM 映射参数
N = 1000;          % 总迭代次数
transient = 200;  % 去除瞬态

% 定义 QO-IM 映射（与 QOIM.m 中一致）
% x_next = 0.2*(-1)^n + a*(phi^2-2)*x;  phi_next = phi + x;
QOIM = @(a, x, phi, n) deal(0.2 * ((-1)^n) + a * ((phi^2) - 2)*x, phi + x);

%% Fig.8(a)：绘制共存吸引子（轨迹）
% 选取两个不同的初始条件
init1 = [0.5, 0.5];
init2 = [0.5, -0.5];

traj1 = zeros(N,2); 
traj2 = zeros(N,2);
traj1(1,:) = init1;
traj2(1,:) = init2;

% 迭代 QO-IM 映射
for n = 1:(N-1)
    [traj1(n+1,1), traj1(n+1,2)] = QOIM(a, traj1(n,1), traj1(n,2), n);
    [traj2(n+1,1), traj2(n+1,2)] = QOIM(a, traj2(n,1), traj2(n,2), n);
end

% 绘制去除瞬态后的轨迹
figure;
hold on;
plot(traj1(transient+1:end,1), traj1(transient+1:end,2), '.r', 'MarkerSize', 2);
plot(traj2(transient+1:end,1), traj2(transient+1:end,2), '.b', 'MarkerSize', 2);
xlabel('x');
ylabel('\phi');
title('Fig.8(a): QO-IM 共存吸引子');
legend('初始条件 [0.5, 0.5]', '初始条件 [0.5, -0.5]');
grid on;
hold off;

%% Fig.8(b)：吸引域绘制
% 构造初始条件网格
nx = 100; ny = 100;
x0_range = linspace(-1, 1, nx);
phi0_range = linspace(-1, 1, ny);
basin = zeros(ny, nx);  % 用数值标记不同类型

maxIter = 500;  % 吸引域判断迭代次数

% 对每个初始条件进行迭代判断
for i = 1:nx
    for j = 1:ny
        x = x0_range(i);
        phi = phi0_range(j);
        diverge = false;
        for n = 1:maxIter
            [x, phi] = QOIM(a, x, phi, n);
            if abs(x) > 1e3 || abs(phi) > 1e3
                diverge = true;
                break;
            end
        end
        if diverge
            basin(j, i) = 3; % 发散：用 3 表示（白色）
        else
            % 对后续 50 次迭代计算 x 的标准差判断是否收敛于固定点
            traj_temp = zeros(50, 1);
            x_temp = x;
            phi_temp = phi;
            for k = 1:50
                [x_temp, phi_temp] = QOIM(a, x_temp, phi_temp, maxIter + k);
                traj_temp(k) = x_temp;
            end
            if std(traj_temp) < 1e-3
                basin(j, i) = 2; % 收敛于点吸引子：用 2 表示（黑色）
            else
                basin(j, i) = 1; % 否则认为进入超混沌吸引子：用 1 表示（青色）
            end
        end
    end
end

% 绘制吸引域（使用自定义色图：1-青色；2-黑色；3-白色）
figure;
cmap = [0 1 1; 0 0 0; 1 1 1];  % 自定义色图
imagesc(x0_range, phi0_range, basin);
set(gca, 'YDir','normal');
colormap(cmap);
colorbar('Ticks',[1,2,3], 'TickLabels',{'超混沌','点吸引子','发散'});
xlabel('x_0');
ylabel('\phi_0');
title('Fig.8(b): QO-IM 吸引域');
